**************************************************************************
*                           VAR.SAS							   *
*                       Joel Hasbrouck							   *
*                        January, 1991							   *
*														   *
*	This file contains the SAS code used to implement the vector	   *
*	autoregressive analysis with random-walk decomposition described	   *
*	in Appendix C of "Using the TORQ Database . . ."				   *
**************************************************************************;

*	Modify the following line to point to your directory;
libname main 'd:\torqprogs\source';
options nocenter linesize=65 mprint symbolgen ps=50;
**************************************************************************
*	Read in Transactions 									   *
**************************************************************************;
data main.trades;
	infile 'd:\torqprogs\sel1ct.asc';
	input symbol $3. tdate 6. ttim 6. ex $1. price 10.
		siz 7. cond $1. g127 3. tseq 6.;
	format ttim time.;	
	run;
proc means fw=6 data=main.trades;
	title 'Stats for Trades';
	run;

**************************************************************************
*	Read in Quotes  	 									   *
**************************************************************************;
data main.quotes;
	infile 'd:\torqprogs\sel1cq.asc';
	input symbol $3. tdate 6. qtim 6. ex $1.
		bid 10. ofr 10. bidsiz 6. ofrsiz 6. mode 5. qseq 7.;
	format qtim time.;
	if ex^='N' then delete;
	run;
proc means fw=6 data=main.quotes;
	title 'Stats for Quotes';
	run;

**************************************************************************
*	Merge Trades and Quotes 									   *
**************************************************************************;
*
*	First lag trades by five seconds;
data main.trades;
	set main.trades;
	ttim = ttim-5;
	run;
data main.varqt;
	merge	main.trades (in=trade rename=(ttim=time))
			main.quotes (in=quote rename=(qtim=time));
	by tdate time;
	retain qm 0;
	if trade then do;
		type = 't';
		if qm=0 then return;
		x0 = 0;
		if price>qm+.0001 then x0=+1;
		if price<qm-.0001 then x0=-1;
		x1 = x0 * siz;
		output;
	end;
	if quote then do;
		type = 'q';
		qm = (bid + ofr)/2.;
		output;
	end;
	format time time.;
	keep tdate time price siz qm x1 x0 type;
	run;

**************************************************************************
*	Next, we have two macro definitions . . .					   *
**************************************************************************;

**************************************************************************
*	The following is a macro to generate SAS code corresponding to 	   *
*	autoregressive lags.  For example:							   *
*	%vargen(gnp,a,0,2)										   *
*	generates the line:										   *
*	a0*gnp + a1*lag1(gnp) + a2*lag2(gnp)						   *
**************************************************************************;
%macro vargen(var,par,lfirst,llast);
%if &lfirst=0 %then
	&par.0*&var;
%else
	&par&lfirst*lag&lfirst(&var);
%if &llast>&lfirst %then %do i=%eval(&lfirst+1) %to &llast;
	+ &par&i*lag&i(&var)
%end;
%mend;

**************************************************************************
*	The following is a macro to estimate a bivariate VAR and		   *
*	perform the random walk decompositions.						   *
*	Parameters:											   *
*	infil	input data set containing r, xs, vr, vxs			   *
*	np		number of autoregressive lags to use				   *
*	nq		number of lags to use in the moving average calculation   *
*	type		"Trades" or "Quotes" (used in annotating the output)	   *
**************************************************************************;
%macro vard(infil,np,nq,type);
*	Compute summary statistics;
proc means fw=6 data=&infil;
	title "&type: Summary Stats (&infil)";
	run;
*	Construct datasets to drive the impulse response computations;
data xsimp rimp;
*	lagged values assumed be zero;
	r = 0;
	xs = 0;
	do i=-&np to -1;
		output rimp;
		output xsimp;
	end;
*	return impulse;
	i = 0;
	r = 1;
	xs = .;
	output rimp;
*	trade impulse;
	i = 0;
	r = .;
	xs = 1;
	output xsimp;
*	future values;
	r=.;
	xs=.;
	do i=1 to &nq;
		output rimp;
		output xsimp;
	end;
*	Estimate the bivariate VAR;
proc model data=&infil;
	title "&type: Estimation of Bivariate VAR (&infil)";
	r = %vargen(r,a,1,&np) + %vargen(xs,b,0,&np);
	xs = %vargen(r,c,1,&np) + %vargen(xs,d,1,&np);
	fit r xs / outest=main.cest outs=main.vcov sur;
	solve / noprint forecast estdata=main.cest data=rimp out=rout;
	solve / noprint forecast estdata=main.cest data=xsimp out=xsout;
	run;
*	Compute the moving average coefficients for the "r" equation;
data rma;
	set rout;
	aa = r;
	set xsout;
	bb = r;
	keep aa bb _lag_;
	output;
data rma;
	array caa(0:&nq) caa0-caa&nq;
	array cbb(0:&nq) cbb0-cbb&nq;
	do i=0 to &nq;
		set rma;
		caa(i) = aa;
		cbb(i) = bb;
	end;
	waa = sum(of caa0-caa&nq);
	wbb = sum(of cbb0-cbb&nq);
	array ssaa(0:&nq) ssaa0-ssaa&nq;
	array ssbb(0:&nq) ssbb0-ssbb&nq;
	sumaa = 0;
	sumbb = 0;
	do i=0 to &nq;
		sumaa = sumaa + caa(i);
		sumbb = sumbb + cbb(i);
		ssaa(i) = sumaa - waa;
		ssbb(i) = sumbb - wbb;
	end;
	do i=0 to &nq;
		aa = caa(i);
		bb = cbb(i);
		if i>0 then waa=0;
		if i>0 then wbb=0;
		saa = ssaa(i);
		sbb = ssbb(i);
		keep i aa bb saa sbb waa wbb;
		output;
	end;

proc print data=rma;
	title "&type: Moving Average Rep of r Equation (&infil)";
	run;
proc print data=main.vcov;
	var r xs;
	title "&type: Disturbance Variance-Covariance Matrix (&infil).";
	run;
**************************************************************************
*	Compute the random-walk and stationary-component variances		   *
**************************************************************************;
data randwalk;
*	read in the variance-covariance matrix;
	do j=1 to 2;
		set main.vcov;
		if j=1 then do;
			vrr = r;
			vrxs = xs;
		end;
		if j=2 then do;
			vxsr = r;
			vxsxs = xs;
		end;
	end;
*	read in the moving average coefficients;
	vs = 0;
	do j=0 to &nq;
		set rma;
		if j=0 then do;
			vw = (waa**2) * vrr + 2 * waa * wbb * vrxs + (wbb**2) * vxsxs;
			vwx = (wbb**2) * vxsxs;
			r2w = vwx/vw;
		end;
		vs = vs + saa**2 * vrr + 2 * saa * sbb * vrxs + sbb**2 * vxsxs;
	end;
	keep vw vwx r2w vs;
	output;
	stop;
proc print data=randwalk;
	title "&type: Results of Random-Walk Decomposition (&infil).";
	run;
%mend;
**************************************************************************
*	End of Macro Definitions									   *
**************************************************************************;

**************************************************************************
*	Perform analysis of Quotes								   *
**************************************************************************;
*
*	First, aggregate over short (5 second) intervals;
data;
	set main.varqt;
	by tdate time;
	retain m 0;
	lasttype = lag(type);
	lasttime = lag(time);
	if first.tdate then do;
		output;
		return;
	end;
	if (lasttype='q' & type='t') or
		(time-lasttime>15)
		then do;
		m = m + 1;
		output;
		return;
	end;
	output;
	drop lasttime lasttype;
data;
	set;
	by m;
	retain x;
	if first.m then x=0;
	if x1^=. then x=x+x1;
	if last.m then output;
data main.varq;
	set;
	qmp = lag(qm);
	if qmp^=. then r = 100*log(qm/qmp);
	x0 = 0;
	if x>0 then x0=+1;
	if x<0 then x0=-1;
	xs = x0*sqrt(abs(x)/100);
	keep r xs;
	run;
%vard(main.varq,5,10,Quotes);

**************************************************************************
*	Perform analysis of Trades								   *
**************************************************************************;
*
*	Pull off trades only;
data;
	set main.varqt;
	if type='t' then output;
data main.vart;
	set;
	pricep = lag(price);
	if pricep^=. then r=100*log(price/pricep);
	xs = x0*sqrt(siz/100);
	keep r xs;
	run;
%vard(main.vart,5,10,Trades);